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Abstract 

We build a minimal, mean-field, model of plasticity of amorphous solids, based upon a 
phenomenology of dissipative events derived, in a preceding paper [A. Lemaitre, C. Caroli, 
larXiv:0705.0823] from extensive molecular simulations. It reduces to the dynamics of an ensemble 
of identical shear transformation zones interacting via the dynamic noise due to the long ranged 
elastic fields induced by zone flips themselves. We find that these ingredients are sufficient to 
generate flip avalanches with a power-law scaling with system size, analogous to that observed 
in molecular simulations. We further show that the scaling properties of avalanches sensitively 
depend on the detailed shape of the noise spectrum. This points out the importance of developing 
a realistic coarse-grained description of elasticity in these systems. 
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I. INTRODUCTION 



Much effort has recently been devoted to the development of theoretical descriptions 
of plasticity of amorphous media. They aim at proposing macroscopic constitutive laws 
consistent with the microscopic information emerging from a wealth of numerical results 
accumulated over the past twenty years. In particular, it is now well established that plastic 
deformation in these systems proceeds via irreversible sudden rearrangements of small clus- 
ters of atoms. Even at zero temperature, these flips occur intermittenly. These empirical 
facts are the basis of two recently proposed phenomenologies: the STZ (shear transforma- 
tion zone) theory [l| and SGR (soft glass rheology) model j^. Both represent the sheared 
disordered solid as a set of spatially random and independent "zones" or "traps" of small 
size embedded in a homogeneous background. These structures are metastable, so that, 
when loaded elastically by the external driving strain, in the absence of noise elastic loading 
would stop at an instability threshold where they flip into an unstressed state. Moreover, 
in both models, structural disorder gives rise to a noise acting in parallel with advective 
loading, thus resulting in intermittent flips occurring before absolute instability is reached. 
An important point is that they introduce noise via Arrhenius factors associated with a 
constant, strain-rate independent effective temperature. 

In an attempt at testing the validity of the assumption concerning flip independence, 
Maloney and Lemaitre ^] later carried out extensive numerical simulations on 2D glasses of 
various sizes L x L in the athermal quasi-static regime (T = 0; vanishing strain rate 7^0), 
hereafter abbreviated as AQS. In this regime, pure elastic, reversible, loading is interrupted 
by randomly spaced discontinuous stress drops associated with the (quasi-instantaneous) 
plastic events. They found that these events can be interpreted as avalanches involving a 
varying number n of elementary zone flips. In the stationary state, the avalanche size n, 
hence the stress drop amplitude Act are broadly distributed, and their averages are system- 
size dependent. The average avalanche size scales roughly as (n) ~ L ['1]. Distributed 
avalanches have also been found by Bailey et al in 3D simulations [sj, with an approximate 
scaling (n) ~ L^/^. 

In the 2D simulations, flips show up clearly as elastic quadrupolar structures in the energy 
and the non-affine displacement field, consistent with their representation as shear trans- 

are known to result from long range 
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interactions. In the present case interzone elastic couplings mediated by the background 
medium are thus likely to be responsible for the avalanche behavior. 

Motivated by this analysis, we have performed in a previous, companion, paper [3] an 
extensive numerical study of the evolution with strain of the non-affine field in a 2D LJ glass. 
It substantiates the importance of interzone elastic couplings and leads to the following 
picture, consistent with that proposed already long ago by Argon et al [8] : via its associated 
quadrupolar field, a zone flip induces, at any other zone site in the system, a shift of the strain 
level whose amplitude and sign depend on the relative position of target and source. These 
shifts may bring some zones past their instability threshold, hence triggering an avalanche. 
For the zones which do not take part in an avalanche, the flip-induced elastic signals act 
as an intrinsic dynamical noise the frequency of which scales as the strain rate 7. The 
statistical study of particle motion shows that this noise dominates largely over fluctuations 
associated with non affine deformations during the elastic episodes separating the plastic 
events. 

In this article we propose a minimal model which incorporates as simply as possible the 
essential features of this phenomenology in terms of a set of identical spatially random zones 
embedded in a 2D elastic continuum, driven by external loading toward their instability 
threshold and coupled via flip-induced quadrupolar elastic fields. 

In Section II we define our model in detail. We discuss the various relevant time scales 
in steady plastic flow, which allows us to identify, for finite systems, a size-dependent quasi- 
static regime where avalanches can be considered instantaneous. We then study numerically, 
in a mean-field approximation for the elastic noise, the steady state dynamics for different 
system sizes. We find that dissipation occurs via broadly distributed flip avalanches, the av- 
erage size of which exhibits a power law scaling with system size: < >~ L^. However, the 
exponent, (3 ~ 0.3, is definitely smaller than the value, of order 1, measured in reference 4|. 

So, in Section III we discuss the various simplifications involved in our modelling of zones 
and of their elastic field. This leads us to test a second version of the model in which we 
assume empirically a gaussian spectrum for the elastic noise. While preserving the existence 
of the avalanche dynamics, this modified model turns out to predict the correct power law 
scaling for < n >. We further show that, in this case, an analytical estimate based on a 
Fokker- Planck- like approximation yields the same prediction for the scaling exponent. 

Although still preliminary, these results lead us to conclude that phenomenologies of 
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FIG. 1: The perturbation due to a localized plastic event corresponds to the elastic response to 
two force dipoles. 

plasticity of amorphous media should definitely incorporate into their basic ingredients the 
strain-rate dependent elastic noise generated by the zone flips themselves. 

II. ELASTICALLY COUPLED ZONES AND AVALANCHES: 
A. The model 

We consider an ensemble of identical zones of size a, shear elastic modulus fi, randomly 
distributed with a fixed density p in a 2D elastic medium of lateral size L. The average 
distance between nearest zones d is fixed and given by : d'^ = L'^/N = a'^/p. The elastic 
state of each zone is characterized by an internal strain e^, which measures the departure 
from its zero stress state. The e^'s lie below a common instability threshold ec- The athermal 
system is driven by external shear at rate e, which advects all e^'s. When a zone reaches ec, it 
disappears while releasing an amount of internal strain Aeo- At the same time another one 
is created, at an uncorrelated position, with zero initial stress (hence zero internal strain). 

During a zone flip the cluster of atoms forming the zone jump into a configuration com- 
patible with the externally imposed strain. This process relaxes the intra-zone stress and, 
at the same time, deforms the surrounding elastic medium. Following Picard et a/ jsl, its 
field can be represented as due to two force dipoles (Figured]). We take for this strain field 
its expression in an infinite medium. At relative position r = (r, 6) from the flipping center: 

Ae(r) = -^^cos4^ (1) 

where 6 is measured from the shearing direction. 

It has quadrupolar symmetry, hence zero average. Aeg and ec can be related using the 
following argument: the amount of stress yuAeo released inside the zone by a flip gives rise to 
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a spatially averaged, macroscopic stress drop Actq = /^fs-AtAeo, where the number k = 0{1) 
depends on the shape of the system. From now on, we assume that k = 1. In stationary 
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state, over a large strain interval, the average number of zone flips is NAe/ec- The associated 
macroscopic plastic stress release NA(Jq must balance the increase of elastic stress, /xAe, due 
to loading. Hence 

^ = = ^ 
ee Na^ a? ^ ' 

The duration tq of a flip is controlled by the time necessary to radiate elastic energy out 
of the zone region, i.e. by radiative acoustic damping, so tq ~ a/cg, with Cg a sound speed. 
This acoustic signal, emitted from site rj, propagates throughout and modifles the strains 
of all other zones j which adjust, over a time ~ tq, to this space-dependent shift Ae(rjj). 
Flip signals thus constitute a noise acting on the e^'s. 

In steady state, the average flip rate in the whole system is 

Tlfiip = ^t'flp = N- = -w- (3) 

If flips occur independently, i.e. in the absence of avalanches, the noise correlation time is tq, 
and the QS regime, where flips can be assumed instantaneous, corresponds to drfup » tq, 
that is to: 

e«e/iip = ^e,— (4) 

For a glass-like system, with e^. ~ 1%, a ~ Inm, zone density a^/(P ~ lO^""^, lateral size 
L ~ 1mm, this yields the loose criterion e << 1. 

Now, a flrst elastic noise signal may drive some e/s beyond ec, hence trigger secondary 
flips, thus initiating an avalanche whose duration Tav is set by sound propagation. For a 
very conservative estimate, we take the average distance between successive flips to be L. 
This leads to a duration Tav =< n > L/cg, with < n > the average avalanche size. It must 
be compared with the average time interval [l^ between avalanches 6tav, given by: 

The QS condition then becomes 

Let us emphasize here an important point, usually ignored in related earthquake models 
centered on the issue of criticallity: the quasi-static range is limited by the acoustic delays 
controlling avalanche spreading, hence shrinks with system size. 
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FIG. 2: n^;, for N = 500. 

We will assume from now on that the QS condition is fulfilled. We thus model the noise 
by instantaneous shifts 5ei of the zone strains. Moreover, we treat our model, now coined "E 
model", in the mean-field approximation, i.e. assume that the 5ej's are independent random 
variables, which amounts to neglecting space correlations between flip centers. That is, 
we take the spectrum n£;(5e) of these noise signals to be that due to a spatially uniform 
distributions of sources truncated at the average distance between nearest zones d. 

ns(5e) = 1^ rdr f dO 6 (Ae(r) - 6e) (7) 



vr (L2 - rf2) Jd 
with Ae(r) given by equation ([1]). 

IlE{Se), plotted on Figure [2l is size-dependent, and has variance 

M,(iV) = -M. (8) 

It exhibits a narrowly peaked structure, associated with distant zones, and vanishes beyond 
the cut-off (2a^Aeo/7r(i2), so that all its moments are finite. Yet, it presents broad, power 
law, tails. For example, its first half moment 

8e, In 



/ d(5e)56nM5e) = ^^ (9) 
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At this stage, the model can be summarized into the following two-step algorithm. Since, 
in the QS regime, the dynamics between instantaneous avalanches reduces to steady drift of 
all ej's at constant loading rate: 



20 



40 
n 



60 



FIG. 3: Distribution vj{n) of avalanche sizes, for our elastic model, and for all system sizes. Insert: 

(n) w{n) vs n/ (n). 

(i) Starting from an initial configuration where all < ec, we first identify eM = Max{ei). 
We shift all ej by Aedrijt = Cc — ^m, and the first avalanche is triggered. Then: 

(ii) Zone M flips, i.e. is removed, while a new one is introduced at zero strain. This first 
flip emits a noise signal which randomly shifts all the other zones: en = Cj + 6ei, where 
the 6eiS are independent and distributed according to IlE{Se). If all en < e^ we are back 
to step (i). Otherwise, an avalanche starts: we count the number qi of zones which flip at 
this stage and are replaced by new, unstrained, ones. 

(Hi) Each zone must now receive zqi signals, which we treat successively. The first of 
these yields a new shifted configuration {£^2} which produces q2 new flips, z is then updated 
to 2; — > 2; — 1 + gi, etc. . . The avalanche stops when z vanishes. Its size is n = 1 + 

B. Numerical results 

In order to study the statistical properties of our model in steady state, we eliminate 
the initial transients (e/cc < 2) and perform ergodic averaging over long strain intervals 
involving ^ 10^ avalanches. 

As shown on Figure El their size distribution WE{n) depends only weakly on the size N 
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for small n( ^ 5). For n ^ 5, it presents a quasi-exponential tail, which broadens noticeably 
with increasing A^. This size dependence reflects into the growth with of the average 
avalanche size (see Figure H]) which we find to fit closely, over the whole iV-range (two 
decades) a power law behavior: 



CtE 



aE = 0.147 



(10) 



When plotting {n)zuE{n) versus n/ {n) (see insert of Figure [3]) it turns out, however, that 
zoE{n) does not obey a simple scaling. Even though power law like decay may be identified 
on a limited, small-n, range, this by no means allows us to conclude to self-criticallity - at 
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variance with a previous claim by Chen et al on a related earthquake model 

We must now confront the above results with those of the molecular simulations. Clearly, 
our highly simplified model shows that long range interzone elastic couplings indeed produce 
broadly distributed avalanches with an average size growing as a power law of system size. 
However, at this stage, agreement is merely qualitative, since (i) The scaling exponent oe 
differs from the simulation value agim — 1/2. (ii) In the molecular simulations |3|, rescaling 
avalanche sizes by their A^-dependent average results in a rather good collapse of data, which 
does not hold in the model. 
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FIG. 4: Average avalanche size (n) vs for the E model. 
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III. DISCUSSION 



A. Approximations of the E model 

The above comparison leads us to discuss in more detail the assumptions involved in our 
minimal model. 

First of all, let us consider more closely our representation of elastic couplings in our 
finite system. When approximating the static elastic propagator by its value for an infinite 
medium, we neglect the finite global elastic recoil which necessarily accompanies, when 
driving at imposed strain, the macroscopic stress release after each flip. This recoil, of order 
should be modelled as a common backward shift, {—^/N) which should be added to 
the noise 6ei. We have rerun the E model for ^ = 1 and 2. We flnd that a flnite recoil leaves 
the power law scaling of {n) (equation ffTOj) ) unchanged, the intuitively expected avalanche 
size reduction showing up only in a slow decrease with ^ of the prefactor. 

Let us now try to list the various simpliflcations underlying our representation of the 
coupled zones. They can be separately into approximations concerning respectively (a) the 
zones themselves and (b) the elasticity of the embedding medium. 

(a) Zones: 

We have taken them to be identical, i.e. to have the same threshold strain ec, and the 
same shear modulus fi, which we have assumed to be constant up to = ec- However, 
clearly, in atomically disordered systems, both fi and ec depend on details of the internal 
structure of the zone and of its immediate vicinity. So, these two parameters are certainly 
distributed about characteristic averages. For example, a signature of the spread of /i is 



the observation in the LJ simulations 



7| of instances in which a zone "overtakes" another 



12| and have checked in ref. 



a 



one during an elastic loading episode. Moreover, we know 
that significant elastic softening occurs near the threshold, where a metastability barrier 
vanishes. 

n 

On the other hand, it was observed in [71] that, frequently, a first fiip does not result in 
the disappearance of the zone. Rather, this persists after a finite strain release, several fiips 
being needed for it to finally "die out". Describing this behavior would demand a multistate 
zone model. 

(b) Elastic couplings: 
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We have represented the embedding elastic medium as a homogeneous continuum, and 
have taken for the elastic field (equation ([1])) its expression for an infinite sytem. We believe 
that this last approximation does not affect avalanche size scaling. Indeed, an elementary 
dimensional analysis for the stain field Ae in a finite L x L box shows that all its moments 
retain the same L dependence as those of the He spectrum (equation ([7j)) used in Section II. 

Note however that the homogeneous continuum approximation itself is probably over- 
schematic, in particular at short distances, where numerical studies of the response to a 
localized force have shown that it is dominated by disorder-induced fiuctuations This 
certainly contributes to a decrease of the tail of noise spectrum. 

Pending more quantitative information about these finer disorder effects, we now propose 
to test the robustness of our minimal model by investigating the avalanche statistics under 
the empirical assumption of a gaussian elastic noise spectrum. 

B. Avalanches in a gaussian noise model 

In this "G model", we choose the gaussian noise spectrum UdSe) to have the same 
variance M2{N) as that (equation ([8])) for the E model. The algorithm is then implemented 
as already described. 

Here again, we find broadly distributed avalanches, whose average size increases with 
The data, shown on Figure El are consistent with the asymptotic (see insert) power-law 
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FIG. 5: Average avalanche size (n) vs N for model G. The dashed line has slope 1/2. 
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FIG. 6: Distribution w{n) of avalanche sizes, for model G and all system sizes. Inserts: {n)w{n) 
vs n/{n). 



Moreover, Figure [6] shows that the rescaling of the size distribution WG{n) as 
(n)^^ f{n/ (n)) leads to good collapse. So, agreement with the results of the molecular 
simulations is much more satisfactory than was the case for the E model. 

We now attempt to clarify how the noise spectrum affects so significantly the avalanche 
behavior. We plot on Figure [7| the strain (ergodic) averages of the distributions pe,g{^) of 
zone strains in steady state for the two models and various system sizes N ranging from 250 
to 32000. For each model, p converges rapidly almost everywhere toward a limit curve. This 
reflects into a very weak size dependence of the macroscopic stress a = 2yu(ej), which varies 
by less than 0.4% (G model) and 0.02% (E model) when N increases from 1000 to 32000. 

While pe and pc are similar in most of the e range, they present significant differences 
in the two regions e ~ and e ~ ec. The peak in p^; results from refeeding zones at e = 
after flips. The larger an avalanche is, the more the corresponding peak is broadened by 
ulterior flips within the avalanche itself. We therefore attribute the washing out of the peak 
for the G model to the fact that it exhibits much larger avalanches. More significant for 
our purpose is the detailed behavior of p near threshold, which reflects avalanche statistics. 



behavior: 





(11) 
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FIG. 7: Steady state distribution p of zone strains vs e/cc for: (a) E model; (b) G model. Inserts: 
blow-ups of the near-threshold region. The arrows indicate increasing values of = 250 x 2'^ 
(g = 0,...,7). 

In particular, = p{€c) is directly related to (n). Indeed, the average flip and avalanche 
rates verify TZ/up = {n)TZav On the other hand, in steady state, the flux of zones which 
cross threhold under the effect of advective elastic loading, i.e. which initiate avalanches, is 
/ = TZav = epcN. So: 

Pc = j\- (12) 

Since the data indicate that (n) diverges for N ^ oo, they also indicate that the absorbing 
boundary condition Pc = should hold asymptotically here. Relation (fT2l) provides a consis- 
tency test of our calculations. We determine pc with the help of a second order polynomial 
extrapolation near ec, with sampling intervals 1.25 x 10~^ec- We find that relation (JT2i) holds 
within 1% for the G model and 3% for the E one. 

In order to try and obtain analytical estimates for avalanche size scalings, we describe 
the evolution of p{e) by the approximate master equation : 

where w is the single flip transition probability 

w{Se) = N-Il{5e) (14) 
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with i the noise distribution, and T = J^^d{Se)w{Se). The delta term, proportional to 
the normalized zone flux f/N = e/ec, accounts for post-flip reinjection and ensures the 
conservation of zone number. 

In this approximation, advection operates between all single flips, which amounts to 
neglecting intra-avalanche time correlations. This leads to an average advective e-shift during 
an avalanche Aeadv ~ av ~ fc(^)/-^! to be compared with the average diffusive broadening 
^^diff ~ \J {n)M2{N), with M2{N) = e^/N the noise variance. Hence, Aeadv/^^diff ~ 



(n) ~ ^ which suggests that our approximation should improve in the large N 

limit. 

Integration of equation (IT^ in steady state yields: 



oo 



^ ^ epc+ I dep{e) I de'w{e ~ e') (15) 



Since w is peaked around zero, in the spirit of the Fokker-Planck approximation, we expand 
p{e) close to ec to flrst order: p{e) ^ + (e — ec)p[- Using pc = {{n)ec)^^, we obtain for the 
avalanche average size: 

-1 



(n) 

where the (semi) -moments 



NMu 
1 + 



(16) 



Mr+ = / d{6e)w{6e)6e' (17) 
Jo 

For both models M2+ = e^/vr^iV, while m[^^ = ec/TrVr^N and Mi^^ is given by equation 
©• 

If p'^ converges towards a flnite value p'^^\ equation ( fT6l) predicts that, for large systems, 

- for the G model: {n) ~ A^^s 

- for the E model: (n) ~ log 

While this prediction accounts satisfactorily for the numerical results for the gaussian 
model, we have checked (see also Figure Hj) that the log scaling is ruled our by our data. 
The reason for this failure is illustrated by the insert of Figure [7] (top). For the E model, 
we flnd that, for increasing A^, p becomes increasingly steep in the near vicinity of ec- In 
the iV-range investigated, we see a marked, non-saturating, increase of | | suggesting a 
possible divergence, higher derivatives increasing even faster. This highly singular behav- 
ior, reminiscent of that analyzed by Chabanol and Hakim [l^ for a related model, clearly 
invalidates the above Fokker-Planck expansion for model E. Conversely, for the G model, 
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we find numerically that v'J^d'^'^ does exhibit convergence, towards ~ —1.3. This regular 
behavior suggests that a Fokker-Planck expansion of the master equation f|T3|) should be 
valid for model G in the large iV limit. We have indeed checked that, for increasing A^, the 
steady state distribution converges rapidly towards the solution of this FP equation. 

So, while even a schematic representation of long range elastic couplings suffices to 
account for avalanches with power law size scaling, this discussion underscores that the 
detailed shape of the noise spectrum is of crucial importance. Indeed, not only does it affect 
the scaling exponent, but, as well, the scaling properties of the distribution of avalanche sizes. 



We consider that the results presented here, though still preliminary, clearly show that 
the dynamical noise due to long range elastic couplings is a key ingredient that must be 
included in phenomenologies of plasticity of amorphous solids. As 7 increases beyond the 
limit of the QS rgime, since avalanches are no longer separable, the spectrum of the flip- 
generated dynamical noise will of course change. The question of its evolution with 7, as 
well as that of its interplay with thermal noise at finite temperature remain for the moment 
completely open ones. Besides, the above discussion indicates two main routes for further 
investigation. 

On the one hand, a more realistic modelization of elastic couplings in the presence of 



structural disorder is needed. Indeed, Leonforte et al |13l . |15[ | have shown that the elastic 
response of amorphous solids self-averages into the continuum elastic response only beyond 
a length scale ^ of order ~ 20 atomic diameters oq. For r < ^, the elastic response is 
dominated by non-affine effects. So, for the existing molecular simulations focussing on 
avalanche dynamics, where L is limited to ^ 50ao, noise tails are very likely to be controlled 
by elastic non-affinity. This issue, as well as that of the statistical distribution of zone 
parameters (such as shear modulus and threshold strain), will demand the development of 
a coarse-grained description of the elasticity of glassy systems. 

On the other hand, our mean-field approximation wipes out from the start the correlation 
anisotropy arising from the quadrupolar symmetry of elementary events, responsible for the 
preferential avalanche orientations observed in the glass simulations of Maloney et al sl 
and Tanguy et al 16|]. In order to evaluate the robustness of the mean-field scalings and 
also start addressing the issue of localization, full simulations of model E in the presence of 
rigid boundaries will be necessary. 
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